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■ ABSTRACT 

, One-dimensional (ID) stellar wind models for hot stars rotating at > 75% of the critical rate show a 

sudden shift to a slow-acceleration mode, implying a slower, denser equatorial outflow that might be 
Cn| ' associated with the dense equatorial regions inferred for B[e] supergiants. Here we analyze the steady 

ID flow equations for a rotating stellar wind based on a "nozzle" analogy for terms that constrain 
the local mass flux. For low rotation, we find the nozzle minimum occurs near the stellar surface, 
allowing a transition to a standard, CAK-type steep-acceleration solution; but for rotations > 75% 
of the critical rate, this inner nozzle minimum exceeds the global minimum, implying near-surface 
supercritical solutions would have an overloaded mass loss rate. In steady, analytic models in which 
the acceleration is assumed to be monotonically positive, this leads the solution to switch to a slow 
acceleration mode. However, time-dependent simulations using a numerical hydrodynamics code show 
that, for rotation rates 75 - 85% of critical, the flow can develop abrupt "kink" transitions from a 
steep acceleration to a decelerating solution. For rotations above 85% of critical, the hydrodynamic 
simulations confirm the slow acceleration, with the lower flow speed implying densities 5-30 times 
higher than the polar (or a nonrotating) wind. Still, when gravity darkening and 2D flow effects are 
accounted for, it seems unlikely that rotationally modified equatorial wind outflows could account for 
the very large densities inferred for the equatorial regions around B[e] supergiants. 

Subject headings: hydrodynamics — stars: early-type — stars: mass loss — stars: winds, outflow — 
stars: emission-line, Be — stars: rotation 
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1. INTRODUCTION sociated with the reduced radiation temperature near 

Hot, luminous, massive stars of spectral types O and B } hc ° quat ° r . can lead to a factor several enhancement 

are generally quite rapid rotators, with inferred surface m the radial mass flux. By itself, this seems made- 

Ctf ■ rotation speeds typically in the range of several 100 km/s, ? uat f to ex P lam equatorial densities estimated to be 

or a substantial fraction of the critical rotation speed (ca. hundreds or even thousands of times the densities^ 

400-600 km/s) at which material at the rotating equato- the polar winds m thesestars (|Zickgraf et al. | [19851; 

rial surface would be in Keplerian orbit. Such luminous i Kraus fc Miroshmchenkol |2Q06|). However, a^eirt_|e- 

stars are also characterized by strong stellar wind out- ries of pape r s by M. Cure and col leagues <|Curi |20fj4j; 

flows, driven by th e line-scattering of the star's contin- I Cure fc R,ial| |2004 fc |Cure et al. | |2005 [) proposes that, for 

uum radiation flux ( Castor et alJ 19751. hereafter CAK). very high near-critical rotation, a switch of the wind out- 

A longstanding question is how such wind outflows are fl ° w to a slo f er ' shallow-acceleration solution can lead to 

affected by the star's rotation, and in particular, whether f imt ^ enhancement m density that, together with the 

this might play a role in the enhanced equatorial density bi-stabihty effect, might reach the equatorial densities 

outflows and/or disks inferred in certain classes of par- inferred in B[eJ supergiants. 

ticularly rapid rotators, e.g., Be and B[e] stars. The present paper aims to understand better the phys- 

For classical Be stars, there is now substantial obser- ical origin of these shallow wind acceleration solutions 

vational evidence (see, e.g., papers in llgnace fe Gavle^l f°, r high rotation rates and to examine critically their 

ESp) that the disks are Keplerian in nature, with very llke jy relevance for explaining dense equatorial disks or 

limited radial outflow; they are thus probably not a di- 0U i> , , . , , . , , 

rect result of feeding by a steady stellar wind (lOwockil Modelli ng rotating hot- storwmds began with the 

2005). However, for the disk and/or enhanced equa- studies by|Friend fc Abbott] (1239, hereafter FA) and by 

torial outflows inferred in supergiant B[e] stars, wind i Pauldrach et & \\ Um hereafter PPK), who extended 

mechanisms seem st ill a viable option. For example , the C ^ K formalism by adding the effect of an outward 

Lamcrs & Pauldrach (1991 and iPelupessv etld] dM ) centrifugal acceleration to one-dimensional (ID) models 

have noted that the "bi-stability" enhancement in opac- 01 ^ ™ nd ° utn ° w ln the equatorial plane Both FA 

ity which occurs for some value of the local surface ef- and PPK independently derived a modified CAK model 

fective temperature (i.e. for B-t ype stars a t ~ 22, 000 ( mCAK ) that relaxes the CAK point-star approx- 

K, see IPauldrach fc Pull p9Ch : iVink et all (HI)) as- lm 1 atlon J ana! Properly accounts for the finite cone-angle 

subtended by the star. Ihey then each found that the 

Electronic address: |toadura@udel.eduowocki@bartol.udel.eduafe ld^ gravity by the Outward Cen- 
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trifugal force tends to increase the mass-loss rate and 
decrease the wind speed. However, for the models com- 
puted, up to about 75% of the critical rotation rate, both 
changes are limited to only a factor of a few and are 
thus insufficient to produce the large equatorial densities 
and low velocities inferred in B[e] supergiants. More- 
over, for still faster rotation, above about 75% criti- 
cal, FA found that the equations for outward accelera- 
tion could no longer be integrated beyond some finite 
radius, and thus they were unable to derive any com- 
plete flow solutions for such rapid rotation speeds. Sub- 
sequent ID models have inves tigate d the role of ma gnetic 
forces dFriend fc MacGregori Il984t iPoe fc Friend! Il986h 
and sound-waves (iKoninx fc HearalT l9921 , but as sum- 
marized bv lBiorkman fc Cassinellil (| 1993) neither mecha- 
nism seems favorable to producing slow, dense equatorial 
outflows. 

More fundamentally, the physical relevance of any such 
ID models may be limited since accounting for latitu- 
dinal flows toward or away from the equator requires 
at least a two-dimensional (2D) treatment. A vivid 
example comes fr om the 2D "Wind Compresse d Disk" 
(WCD) model of iBjorkman fc CassinellH (fl993h . which 
argues that conservation of angular momentum should 
tend to channel material from higher latitudes toward 
an equatorially compressed disk flow. If one assumes a 
purel y radial driving forc e, 2D hydrodynamical simula- 
tions (jOwocki et al.1119941 ) confirm the basic WCD effect, 
but show that, depending on whether material reaches 
the equator above or below some "stagnation point" , it 
either drifts outward or falls back toward the star. Such 
simultaneous infall-plus-outflow behavior is not possible 
in a steady ID model, but is a perfe ctly natural outcome 
in a 2 D simulation. Furthermore, ICranmer fc Owockil 
(1995) showed that, when computed from proper angle 
integration of intensity from the rotationally distorted 
stellar surface, the line- force also has nonz ero, nonradial 
compo nents in both azimuth and latit ude. lOwocki et al.l 
(119961 see also lPetrenz fc Puis! ((2000)) showed that such 
nonradial line forces can inhibit the formation of a WCD. 
Finally, equatorial gravity darkening can actually reduce 
the wind mass flux from the equator, and so lead to an 
equatorial wind density that is lower, not higher, than 
near the poles. 

Despite this likely i mportance of such 2D effects , 
several recent a nalyses (|Curel l2004t ICure fc Rial 12004 
ICure et aLll2005f ) have reexamined the ID equatorial flow 
models of FA, with a particular focus on the failure to 
obtain monotonically accelerating wind solutions for ro- 
tation above about 75% of the critical speed. In partic- 
ular, these papers argue that for such very high, near- 
critical rotation, the wind solution can switch to an al- 
ternative mode, characterized by a much slower out- 
ward acceleration. Together with a moderately enhanced 
mass flux, the resulting lower-speed outflow then im- 
plies a substantial enhancement in density, relative to the 
standard CAK, steep acceleration applicable at higher 
latitudes. When combined with parameterizations in- 
tended to mimic a "bi-stability" enhancement in the 
line-d r iving (lLamers fc Pauldrach|[l991t iPelupessv et al.1 
12000ft , ICure et al.1 (|2005f ) predict equator-to-pole density 
contrasts of the order of 10 2 — 10 4 . 

Assessing the physical relevance of such claims for 
understanding B[c] stars will eventually require proper 



account of the multi-dimensional effects noted above. 
Nonetheless, to provide a solid basis for such multi- 
dimensional models, it is important to have a clearer, 
dynamical understanding of these novel ID slow outflow 
solutions. Combining analytic studies with numerical hy- 
drodynamic simulations, this paper examines the reality 
and physical origin of the shallow wind acceleration so- 
lutions for high rotation rates, with emphasis on their 
possible relevance to disk outflows from B[e] supergiants. 

We begin (§[2]) with a basic review of the general time- 
dependent wind equations, together with their CAK- 
type, steady-state solutions in a non-rotating wind. To 
provide a physical basis for extending these steady mod- 
els to include rota tion, we first (§151) apply a s imple "noz- 
zle" analysis (see lHolzerlll977t lAbbottl 11980ft . originally 
developed to study winds driv en from luminous accre- 
tion disks ((Perevra et al.1 120041 ). With a few judicious, 
yet quite reasonable approximations (e.g., neglecting gas 
pressure terms by taking the zero-sound-speed limit; us- 
ing a beta velocity law to evaluate the finite-disk correc- 
tion as an explicit function of radius), it is possible to 
obtain simple integrations of the equation of motion and 
to study the scalings of the mass loss rate with rotation, 
as well as the switch from steep to shallow acceleration 
solutions beyond a threshold rotation rate. To test the 
validity of these simple analytic solutions, we next apply 
a numerical hydrodynamics code (§|4|) to evolve ID ro- 
tating wind models to asymptotic steady states f§ !5.1[) . 
Results confirm a transition to slower acceleration at very 
high rotation (above about 85% of critical), but also show 
a new class of non-monotonic "kink" solutions that apply 
for moderately fast rotation (ca. 75-85% of critical) . We 
then examine (§ !5.2ll5~3"l) the time evolution of solutions in 
various rotation domains, with emphasis on the kink so- 
lutions, and on a peculiar transition case (85% of critical 
rotation), characterized by an initial wind overloading 
followed by a flow stagnation and eventual reaccretion 
of material onto the star. We conclude with a summary 
and outlook for future work (§[6]). 

2. GENERAL FORMALISM FOR LINE-DRIVEN MASS-LOSS 

2.1. ID Time- Dependent Equations of Motion 

In this paper, we examine ID radiatively driven out- 
flow in the equatorial plane of a rotating star. For the 
general time-dependent simulations discussed in §3, the 
relevant equations for conservation of mass and radial 
component of momentum have the form 

dp 1 d(r 2 pv) 
_ + 

I 111' f T W (I — I I 

■ Alines, (2) 







dv dv 
dt dr 



2 Q r 

18P GM*(l - T e 
p dr r 2 



(1) 



where r and t are the radius and time, and p and v are the 
mass density and radial component of the velocity. The 
body forces here include the outward radiative accelera- 
tion from line-scattering, <?ii nC s, and an effective inward 
gravitational acceleration GM,(1 — r e )/r 2 , reduced by 
the outward continuum radiative force from scattering 
by free electrons, as accounted for by the Eddington pa- 
rameter r e = n e L*l 'AttGM^c. For the centrifugal term, 
v^/r, we avoid explicit treatment of an azimuthal mo- 
mentum equation by assuming simple angular momen- 
tum conservation (which is a good approximation in the 
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supersonic flow domain considered here), yielding then 
for the azimuthal speed 

V<t> = l>rot — j (3) 

where v mt is the rotation speed at the star's equatorial 
surface radius i?* . For simplicity, we also avoid a detailed 
treatment of the wind energy balance by assuming an 
isothermal outflow, for which the pressure is written as 
P = pa 2 , where a is the (constant, isothermal) sound 
speed. 

2.2. Steady-State Equations with Rotation 

For the simplified case of a steady state, the time- 
dependent terms vanish (d/dt — 0), yielding for the 
steady acceleration 



dv 
dr 



GM«(l-r e 



v 2 B 2 



■fflix, 



(4) 

p dr 



The steady form for mass conservation implies a con- 
stancy for the overall mass loss rate, M = 4irpvr 2 . Using 
this to eliminate the density p, the equation of motion 
([3]) takes the form 



v 



dv 

! 

dr 



GM.(l-r e ) | vl 



2a 2 



-fflix. 



(5) 

The square-bracket factor on the left-hand-side allows for 
a smooth mapping of the wind base onto a hydrostatic 
atmosphere below the sonic point, where v < a. How- 
ever, in radiatively driven winds the pressure terms on 
the right-hand-side are generally negligible since, com- 
pared to the gravitational acceleration term that must 
be overcome to drive a wind, these are of order w s = 
(a/v csc ) 2 w 0.001, where v csc = y/2GM*(l - T e )/R* is 
the effective escape speed from the stellar surface radius 

Since the key to a stellar wind is to overcome gravity, it 
is convenient to define a dimensionless equation of motion 
that scales all accelerations by the effective gravity, 



(1 - W s /w)w' = -1 +L0 2 (l — x) + ru nes + 



4w» 



(1-x) 



(6) 



3iincs?' 2 /GM»(l — r e ), and the gravitation- 
ally scaled inertia is w' = dw/dx = r 2 v(dv/dr) /GM t (l- 
T e ). The independent variable here is the inverse radius 
coordinate x = 1 — R*/r, while the dependent variable is 
the ratio of the radial kinetic energy to effective surface 
escape energy, w = v 2 /v 2 sc . Gas pressure effects are ac- 
counted for by terms containing w s = a 2 /v 2 sc , while cen- 
trifugal effects from rotation are characterized in terms 
of the ratio of the equatorial rotation speed to critical 
speed, u> = u r ot/w C rit = V% v I0 t / v eS c, under the assump- 
tion that the wind material conserves its surface value of 
specific angular momentum, rv${r) = v ro tR*. 

Within the CAK formalism for driving by scattering of 
a point-source of radiation by an ensemble of lines, the 
deshadowing of optically thick lines by the Doppler shift 
associated with the wind acceleration gives the scaled 
radiative acceleration ri; ncs a power-law dependence on 
the flow acceleration w'. 



(7) 



where a is the CAK power index. Here we have elimi- 
nated an inverse dependence on density p in favor of the 
mass loss rate M — 4irr 2 pv, with the line force constant 
thus defined by 



C = 



1 



1 



Mc 2 



(8) 



.i-r e 

with L* the stellar luminosity, and r e the Eddington 
parameter for the gravitationally scaled radiative ac- 
celeration from electron scattering opac ity, K e (cm 2 /g) 
(lLamers fe Cassinellilll999t IOwockill200l . We have also 
used the lGavlevi (|1995 ) Q notation for the overall normal- 
ization of the line opacity. Note that, for fixed sets of stel- 
lar parameters {L Sfl M Sr ,T e ) and line-opacity (a,Q), the 
constant C scales with the mass loss rate as C oc 1/M a . 

As already noted, the smallness of the dimension- 
less sound-speed parameter w s implies that gas pressure 
plays little role in the dynamics of any line driven stellar 
wind. Hence, to a good approximation, we can obtain 
accurate solutions by analyzing the much simpler limit 
of vanishing sound speed a oc yjwl — * 0, for which the 
line-driven- wind equation of motion reduces to 

w' = -l+uj 2 {l-x) + Cw' a . (9) 

2.3. Classical CAK Solution for a Point-Star 

Let us first review the standard CAK solution without 
rotation, setting ui = 0. Note then that since the pa- 
rameters T e and C are spatially constant, the solution is 
independent of radius. For high M and small C there 
are no solutions, while for small M and high C there are 
two solutions. The CAK critical solution (denoted by 
the subscript c) corresponds to a maximal mass loss rate, 
which requires a tangential intersection between the line- 
force Cw' a and the combined inertia plus gravity 1 +it/, 
for which 

aC'cw'^ 1 = 1 , (10) 
and thus, together with the equation of motion ([9|), we 
have 

rv 

(11) 



with 



CI 



1 - a 



1 



(12) 



a«(l-ar " 

Using eq. ([5]) , this then yields the standard CAK scaling 
for the mass loss rate, 



CAK — 



(l-a)/a 



(13) 



Moreover, since the scaled equation of motion ([9]) has no 
explicit spatial dependence, the scaled critical accelera- 
tion w' c applies throughout the wind. This can therefore 
be trivially integrated to yield 

w(x) = w(l)x , (14) 

where w(l) — a/ (I — a) is the terminal value of the scaled 
flow energy. In terms of dimensional quantities, this rep- 
resents a specific case of the general "beta" -velocity-law, 



v(r) 



i?* 



(15) 



where in this case (3 = 1/2, and the wind terminal speed 
scales with the effective escape speed from the stellar 
surface, Voo = v csc ^a/{\ - a). 
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wind. For a rotating star and wind, the fdcf can be- 
come even more complicated, modified by the rotational 
shear of the wind outflow, and by the oblateness of the 
star, and possibly also by th e equatorial gravity dark- 
ening of the source rad iation (jCranmer fc Owockil Il995l : 
iGavlev fc Owockil 12000). However, for simplicity, let us 
nonetheless base our analysis on the spatially explicit 
form obtained by assuming a canonical (3 = 1 velocity 
law (p~5|) within the finite disk factor p^|) for a simple 
spherical expansion. In the zero sound-speed limit, the 
scaled equation of motion ([9]) can now be written in the 
form, 



Fig. 1. — Spatial variation of the finite-disk-correction factor /, 
plotted vs. scaled inverse radius x = 1 — R t /r, for CAK exponent 

= l/2 and various velocity-law exponents, /3 = 0.5, 1, 2, and 
5. The horizontal dashed curves denote the unit correction that 
applies at the point of isotropic expansion (where a = dlnv/dlnr— 

1 = 0), and at large distances where the star approaches the point- 
source form assumed in the original CAK model. 



2.4. Finite-Disk Form for the CAK Line-Force 

The above analysis is based on the idealization of radi- 
ally streaming radiation, as if the star were a point source 
at the origin. This was the basis of the original CAK 
wind solutions, although they did already identify (but 
did not implement) the appropriate "finite-disk correc- 
tion factor" (fdcf) to account for the full angular extent 
of the star (see CAK eq. [50 ]), 



-l+LU 2 (l~x)+fC c 



(17) 



m = 



(i 



\l+a 



2\l+a 



(l + a)*(l + *) a (l-iJ*) 



(16) 



where we have normalized the line force in terms related 
to the "point-star" CAK model, with m = M/Mqak, 
the ratio of the mass-loss rate to the point-star CAK 
value. Note then that for the non-rotating (ui = 0), 
point-star (/ = 1) case of the classical CAK model, the 
critical solution (with maximal mass loss) is given by 
m = 1 and w(x) — a/ (1 — a)x. As noted in i)2.3l this 
implies a CAK mass loss rate M — McAK an d a ve- 
locity law, v(r) = VooyJ\ — R*/r, with terminal speed 

Woo = y/a/(l- a)v csc . 

To analyze models with rotation, a particularly con- 
venient case is to take a — 1/2, for which the equation 
of motion (fT7|) (using eq. (fT2|) for C c ) becomes a simple 
quadratic in 



'/ - 2fy/w'/m + g(x) = , 



(18) 



with u* = \/T — R 2 Ir 2 the cosine of the finite-cone an- i. c ■ i_ j c j a a- u 

i r ;i ; ii i* i i 7i ; n -, ,,T, where for convenience, we have defined a rotationally re- 

rein <-v4- 4- h n n+nl I rlinlr nnrl i-r W In /ii / W Iti m I \A/hnn ' _ *■' 

duced gravity as g(x) = 1 — or (1 — x). We can then solve 
for a shallow (-) and steep (+) acceleration solution, 



gle of the stellar disk, and a = dh\v / d\nr — 1. When 
this factor is included to modify the point-star CAK line 
acceleration (from eq. [7]), its complex dependence on 
radius, velocity, and velocity gradient complicates the 
solution of the full equation of motion. Full solutions 
derived independently by FA and PPK yield a some- 
what reduced mass-loss rate Mfd ~ Mcak/(1 + a) 1 ^" 
and higher terminal speed Uqo f=a 3u esc - 

But if we approximate the wind velocity law by the 
simple "beta-law" form of eq. (fl~5|) , then the fdcf can 
be evaluated as an explicit spatial function. Figure [T] 
illustrates the resulting variation of / with the scaled co- 
ordinate x for a = 1/2 and various values of (3. Note that 
the overall form is quite similar for all cases, increasing 
from a surface value /* = f(R*) = l/(l+a) to past unity 
at the isotropic expansion radius (where dv/dr — v/r), 
r/R* = (1+0) [corresponding to x = (3/(1 + /?)], and 
eventually returning asymptotically to unity from above 
at large radii (x — > 1). 

In the steady wind analysis in the next section, we thus 
choose the canonical value (3 = 1 to represent the fdcf as 
an explicit spatial function (PPK). 

3. STEADY-STATE SOLUTIONS FOR ID MODELS OF 
ROTATING, LINE-DRIVEN STELLAR WINDS 

3.1. Nozzle Analysis for Steady Wind Acceleration 

Let us now examine how the combined effects of the 
fdcf and rotation alter the classical CAK result. Note 
that we are ignoring here gravity darkening and oblate- 
ness effects, as well as any "bi-stability" in the line- 
driving parameters between the polar and equatorial 



w' ± (x) 



g(x)n(x) 



1 ± y/l — m/n(x) 



with the "nozzle function" , 



n(x) 



fix) 



m 2 

g(x) 1 — uj 2 (1 — x) 



(19) 



(20) 



The significance of this nozzle function stems from its 
appearance with the mass loss rate m within the sq uare- 
root discriminant (cf. Laval nozzle; lAbbottJ (|1980D 'I. In 
particular, we can readily see that maintaining a numer- 
ically real flow acceleration requires 1 a mass loss rate 
m < min[n(x)]. As such, the location of the global mini- 
mum of this function (the smallest nozzle "throat" ) rep- 
resents the critical point that sets the maximal allowed 
value of the mass loss rate, m = min[n(a;)], that is con- 
sistent with a monotonically accelerating outflow. 

Figure [5] plots n(x) vs. x for various rotation rates 
ui, using a (3 = 1 velocity law to obtain a spatially ex- 
plicit approximation to the fdcf. Note that for no or 
low rotation (about lu < 0.75), the minimum of the noz- 
zle function is less than unity, and occurs at the stellar 

1 Actually, this restriction really stems from our CAK scaling of 
the line-force with w la (in this case VuF), which requires a strictly 
positive acceleration, w > 0. But if we provide a backup scaling 
for negative accelerations, then "overloaded" situations, for which 
the square-root discriminant in eq. H19H beco mes negative, simply 
lead to an abrupt switch, a so-called "kink" (Cranmcr & Owockil 
1996), to a decelerating solution. See i[4l and l6l for further details. 
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past the threshold rate atw~ 0.75, the solid curves show 
an abrupt shift from steep acceleration to shallow accel- 
eration, with the mass loss saturating to the point-star 
CAK value, m = 1. The dashed curves show extrapo- 
lated results if the local nozzle minimum at the surface 
is instead used to set flow conditions; the mass loss in 
this case is set by the scaling mo in eq. (|2ip. and the 
terminal speed is derived by assuming a pure gravita- 
tional coasting for all radii with n{x) < mo- The data 
points again compare corresponding results for the full 
dynamical simulation, as described further in the next 
sections. 



Fig. 2. — Nozzle function n(x) plotted vs. scaled inverse radius 
x = 1 — R*/r, for various rotation rates w = 0, 0.5, 0.6, 0.7, 0.8, 0.9, 
0.95, ranging from lowermost to uppermost. All curves use a f3 = 1 
velocity law in evaluating the fdcf. The horizontal dashed line at 
unit value represents the nozzle function for the CAK point-star 
model, with n = / = g = 1. 



surface, x — 0. This allows the flow to transition to 
a super-critical outflow directly from the static surface 
boundary condition w(0) = 0, following the steeper, plus 
(+) root for the acceleration in equation (TT9"]) . but with 
a mass loss rate less than the point-star CAK value, 



m = m = n(x = 0) 



f 2 (x = 0) 
1 -w 2 



4/9 



(21) 



Note that the factor 4/9 in the numerator is just the 
a = 1/2 value for the zero-rotation, finite-disk-corrected 
mass loss scaling derived by FA and PPK, 



m fd 



Mfd 

MCAK 



- fl /a = 



1 



(l + a)V« 



(22) 



By contrast, for large rotation rates (about cu > 0.75), 
this nozzle minimum is unity, and occurs at large radii, 
x = 1; satisfying the static surface boundary condition 
now implies that the flow at all finite radii should remain 
sw&-critical, following the shallower, minus (-) root for 
the acceleration in equation (fl9|) . now with a mass loss 
rate just equal to the point-star CAK value, m = 1. 

This thus provides the basic explanation for the 
switch from steep to shallow accelerations inferred by 



switcn irom steep 
iCure et all (I2005I ). 



3.2. Nozzle Solutions for the Velocity Law in ID 
Rotating Winds 

The associated wind velocity laws can be obtained 
by simple numerical integration of equation (|19p from 
a static boundary w(Q) = 0, following either the steep 
or shallow solution, depending on whether the rotation 
rate is high enough to shift the critical point (where n{x) 
has its absolute minimum) from the surface (x = 0) to 
large radii (x — 1). Figure [3^ plots the resulting velocity 
laws for selected slow vs. rapid rotation rates, yielding 
respectively the steep vs. shallow types of flow solution. 
The dashed curve in Figured plots the escape speed v esc 
as a function of x, showing that these winds are capable 
of escaping the star. The right panel compares results 
from full dynamical simulations described below. 

Figure |4] illustrates the associated terminal speed and 
mass loss rates for these solutions (solid lines), plotted as 
a function of rotation rate u>. As the rotation increases 



4. SPECIFICATIONS FOR NUMERICAL HYDRODYNAMICS 
SIMULATIONS 

The above nozzle analysis provides a helpful frame- 
work for understanding the nature of flow solutions from 
a rotating wind model, but is based on some key simpli- 
fications, e.g. neglect of gas pressure (inclusion of which 
would lead to a wind that is not super-critical directly 
from the static surface boundary) and an approximate, 
spatially explicit form for the finite-disk correction factor. 
Moreover, it implicitly assumes that the derived solutions 
are the only relevant stable, attracting steady states for 
the rotating wind outflow. To test the validity of these 
simplifications and assumptions, let us now examine the 
time evolution of analogous ID flow models, including 
also both the finite gas pressure and the dynamically 
computed finite disk correction factor. Our specific ap- 
proach here is to use a numerical hydrodynamics code 
to evolve a ID time-dependent model of the equatorial 
plane a line-driven stellar wind from a rotating star to- 
ward an equilibrium steady-state for the resulting flow. 
The results of these simulations can then be compared to 
those predicted by the steady-state nozzle analysis pre- 
sented above. 

The numerical models presented here were com- 
puted using a piecewise parabolic method (PPM; 
ICollela fe~ Woodward 1984) hydrodynamics code called 
VH-1, originally developed at the University of Vir- 
ginia (J. Blondin). The basic code was modified for the 
present study to include radiative driving terms, solv- 
ing the time-dependent equations for ID spherical out- 
flow (TT]) and ([SJ> - The spatial mesh uses n r = 600 ra- 
dial zones extending from the base at the stellar surface, 
r\ = R ni i n = R*, to a maximum r 60 o = i?max = 100i?*, 
with the zone spacing starting at Ari = 6.18 x 10~ 5 i?* 
and then increasing by 2.5% per zone out to r35o = 15i?*, 
after which it remains constant at Ar = 0.333 i?* to the 
outer boundary. Tests with double the grid resolution for 
selected cases give similar results to the above standard 
values. 

The parabolic method requires flow variables to be 
specified in a phantom zone beyond each boundary. At 
the outer radius, we assume supercritical outflow, with 
boundary conditions set by simple flow extrapolation as- 
suming constant gradients. This is justified because, 
as discussed in the Appendix, when finite sound-speed 
terms are included, the critical point for even the shallow 
acceleration solutions should be well within our assumed 
outer boundary radius of i? max = 100i?„. 

At the inner boundary, the velocity in the two radial 
zones below i m i n is set by constant-slope extrapolation, 
thus allowing the base velocity to adjust to whatever is 
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Fig. 3. — (a) Flow speed over escape speed, v/v eBC = y/w, plotted vs. scaled inverse radius x = 1 — R*/r, as derived from a nozzle 
analysis, using a steep acceleration for no or modest rotation, oj = 0, 0.5, 0.6, 0.7, 0.74 and shallow acceleration for rapid, near-critical 
rotation, ui = 0.75, 0.8, 0.85, 0.9, 0.95. The dashed curve shows the escape speed as a function of x. b) Same plot as in (a), but as found 
using asymptotic states of full hydrodynamical simulations. Note again the steep, supercritical accelerations for no or modest rotation, 
u) = 0, 0.5, 0.7, and shallow, subcritical accelerations for rapid, near-critical rotation, w =0.86, 0.9, 0.95. However, note also the "kink" 
solutions present for u) = 0.75, 0.8, 0.82, 0.84 and the collapsed solution for w = 0.85. 
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Fig. 4. — Upper solid curve: Terminal flow speed over escape 
speed, foo/fese = \Au(l), plotted vs. rotation rate uj, showing 
the shift from fast to slow wind as rotation rate is increased past 
ca. uj = 0.75. Lower solid curve: Mass loss rate in units of point- 
star CAK value, again plotted vs. rotation rate uj, showing the 
saturation at the CAK mass loss rate for rapid rotation, ui > 0.75. 
Dashed curves show the continued steady decrease in Veo/v eBC and 
increase in fn if the local nozzl e mi nimum at the stellar surface 
is used to set flow conditions (5 13.20 . The circles (i>oo/fesc) and 
triangles (m) show corresponding results from full hydrodynamical 
simulations. 



appropriate for the overlying flow (|Owocki et al.lll994l) . 
This usually corresponds to a subsonic wind outflow, al- 
though inflow at up to the sound speed is also allowed. 
The base density is fixed at p — 8.709 x I0~ 13 g cm~ 3 , 
a value chosen because, for the characteristic wind mass 
fluxes of these models, it yields a steady base outflow 
that is moderately subsonic. A lower-boundary density 
much smaller than this produces a base outflow that is 
supersonic, and thus is unable to adjust properly to the 
mass flux appropriate to the overlying line-driven wind. 
On the other hand, a much larger base density makes the 
lower boundary too "stiff" , leading to persist ent oscilla- 
tions in the base velocity (jOwocki et al.lll99"4T) . 

These time-dependent simulations also require setting 
an initial condition for the density and velocity over the 



entire spacial mesh at some starting time t — 0. For this 
we generally use a standard, finite-disk-corrected CAK 
wind, computed by relaxing a ID, non-rotating simu- 
lation to a steady state; however, for selected models 
with moderately rapid rotation, we also explore using 
a slow-acceleration initial condition (see § 15. 3[) . From 
the assumed initial condition, the models are advanced 
forward in time steps set to a fixed fraction 0.25 of the 
Courant time. 

Our version of VH-1 is set up to operate in CGS units, 
requiring specification of physical values for the basic 
parameters for both the star (e.g., mass, radius, lumi- 
nosity) and wind (e.g., CAK k, a, and S). Building 
upon our earlier studies of Be stars, the specific pa- 
rameters chosen here are for a main sequence B star 
with mass M* = 7.5M Q , radius i?* = 4i? Q , luminosity 
L* = 2310L Q , and temperature T = 2 x fO 4 K; but we 
have also explored models with parameters appropriate 
for supergiant B[e] stars. These stellar parameters imply 
an Eddington parameter of T e — 0.008, an isothermal 
sound speed of a — f6.6 km/s, and escape and critical 
speeds of v csc = 845 km/s and u cr it = 597 km/s. We also 
assume a CAK power-law index of a = 1/2 and cumu- 
lative line-strength parameter Q — f533 (|Gavlevi rf 995). 
The value of S has been set to zero in all simulations. 

In any case, for a given choice of the CAK power-law 
index a, we find that results are largely independent of 
the specific physical parameters when cast in appropri- 
ately scaled units, normalizing for example radius by the 
stellar radius i?*, velocity by the wind terminal velocity 
Woo (which in turn scales with the stellar surface escape 
speed), time by the characteristic flow time i?*/«oo, and 
mass loss rate in terms of the classical (point-star) CAK 
value given in eqn. (|I3|) . To faciliate comparison with 
the analytic nozzle analysis in §§2-3, we again choose 
a = 1/2, and plot all simulation results using the above 
scalings. 

We further note that essentially all the key VH- 
I results reported here were very well reproduced 
by a completely independent, simple dimensionless 
hydrodynamics code developed by one of us (AF; 
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iFeldmeier fc Nikuttal[200l . 

Finally, in such time-dependent simulations of line- 
driven winds, one must also supply a generalized scal- 
ing for the line-force that applies in the case of non- 
monotonic flow acceleration. In general this requires tak- 
ing into a ccount non-local couplings of the line-transfer 
(see, e.g.. IFeldmeier fc Nikuttal 12(306). but we wish here 
to retain the substantial advantages of using a purely 
local form for the line-driving. Noting that a nega- 
tive velocity gradient implies a prior line resonance that 
shadows radial photons from the star, a lower limit 
would be just to set gu nc — whenever dv/dr < 0. 
On the other hand, since forward scattering can sub- 
stantially weaken any such shadowing by a prior reso- 
nance, an upper limit would be to compute the local line 
force using the absolute value of the velocity gradient, 
giino <x \dv/dr\ a . As a simple compromise between these 
two extremes, we choose here a scaling that truncates 
the radial velocity gradient to zero whenever it is nega- 
tive, i.e. dv/dr — > max(rfu/rfr, 0) . For a point-star model 
with radially streaming radiation, this would give a zero 
line- force (since dv/dr < 0), but when one accounts for 
the lateral expansion v/r within the finite-disk correc- 
tion factor, it leads to a line force in which the usual 
dependence on radial velocity gradient is replaced by a 
dependence on the expansion gradient, 

Siinc cx (v/r) a . (23) 
This leads to a line acceleration that is intermediate be- 
tween the underestimate and overestimate of the two 
more extreme scalings. 

5. COMPARISON WITH TIME-DEPENDENT 
HYDRODYNAMICAL SIMULATIONS 

5.1. Asymptotic Steady- States of Time- Dependent 
Simulations 

In our basic parameter study, each simulation is run 
using the same initial condition, spatial mesh, bound- 
ary conditions, base density, etc., with the only varia- 
tion being the rotation rate uj, which is set to specific 
values ranging from 0.1 to 0.97. With only one excep- 
tion (for the uj = 0.85 case, which turns to be a rather 
pathological value; see ? 15. 2[) . all simulations asymptoti- 
cally relax to a well-defined steady state. Moreover, for 
both moderate rotation (uj < 0.75) and high rotation 
(uj > 0.85), these asymptotic states agree remarkably 
well with the predictions of the above nozzle analysis. 
Figure [3)3 shows the velocity laws for these final states, 
scaled in the same form used in figure 3a for the nozzle- 
analysis results. Note that both figures 3a and 3b show a 
steep acceleration for no or modest rotation (uj = 0-0.7), 
and shallow acceleration for rapid, near-critical rotation 
(uj = 0.86, 0.9 and 0.95). 

However, for the moderately high rotation-rate cases 
uj = 0.75, 0.8, 0.82 and 0.84, note also the appearance of 
a new class of kink solutions, characterized by an abrupt 
shift to a decelerating or "coasting" flow beyond a well- 
defined "kink radius" rtink- These kink solutions thus 
represent a kind of intermediary final state of the time- 
dependent simulations in the parameter ranges 0.75 < 
uj < 0.85, effectively smoothing the abrupt jump from 
fast to slow acceleration solutions expected from steady- 
state analyses. The formation of such kinks, and their 
underlying physical cause, are discussed further in the 
section below (§ I5.2[) on time evolution. 



Figure [4] shows that the fully dynamical results for 
the scaled ratio of right boundary speed (circles) and 
CAK scaled mass loss rate (triangles) are generally in 
good agreement with the predictions of the simple nozzle 
analysis (solid curves), with the modest, ca. 10% differ- 
ences likely attributable to inclusio n in the simulations 
of a s mall but finite sound speed (|Owocki and ud-Doulal 
12004) . However, for rotation rates 0.75 < u> < 0.85, the 
dynamical results tend to follow the dashed curves of 
the extended nozzle analysis, representing extended fast 
solutions rather than the abrupt shift to slow solutions 
indicated by the solid curves. As discussed above, this 
range of rotation rates is characterized by kink solutions. 
To understand better this development of fast vs. slow 
vs. kink solutions, let us now examine the time evolution 
of the simulations toward asymptotic states. 

5.2. Time Relaxation of ID Rotation Models 

For the specific rotation rates w = 0.7, 0.8, 0.84, & 
0.90, which span the parameter range between fast and 
slow acceleration solutions, Figure[5]uses gray-scale plots 
of the mass-loss rate (in units of the point-star CAK 
value) to illustrate the time relaxation from the CAK 
initial condition (set to a non-rotating, finite-disk cor- 
rected CAK model) to an asymptotic steady-state. Time 
is given in units of the flow time t — i?*/foo, and the final 
values of the scaled mass-loss rate are indicated above the 
associated plot. Note that there are distinct differences 
in the time evolution of each model, with the more rapid 
rotation cases characterized by a longer relaxation time, 
but the temporal and spatial constancy of the final mass 
loss rates illustrate the steady nature of the asymptotic 
solutions. 

Despite this temporal and radial constancy in the mass 
flux, the asymptotic states of the models with uj = 0.8 
and 0.84 include a kink, or abrupt discontinuity in their 
velocity gradient. To illustrate the formation of these 
kinks, the gray-scales in figure [5] show the evolution of 
the scaled velocity gradient. In both cases, the abrupt 
switch from a steep acceleration solution to a decelerating 
solution is apparent from the sharp transition from pos- 
itive to negative velocity gradient, with the kink radius, 
Hank, indicated above each plot. Since, despite the radial 
discontinuity, the velocity gradient contours all become 
constant in time, we see clearly that the kink solutions 
are indeed perfectly valid steady-states. 

Note, however, that the kink radius shifts from rkmk = 
4.18i?* for the u = 0.8 model to r kin k = 2.12i?* for the 
more rapid-rotation, uj = 0.84 model. From the velocity 
law plots in figure [3Jd, this also implies that the more 
rapid rotating case has a more extended deceleration re- 
gion, and thus ends up with a much lower final speed. We 
can thus anticipate that a somewhat faster rotation rate 
should give an even lower kink radius, with the more ex- 
tended decelerating region leading to an even lower final 
speed, or perhaps even to a flow stagnation (zero veloc- 
ity) at a finite radius. In particular, note from figure [4] 
that the dashed extrapolation curve suggests that the on- 
set of such flow stagnation should occur near a rotation 
rate of uj ps 0.85. 

In fact, our numerical simulations do show that models 
near this rate have a quite pathological behavior. This 
is illustrated in figure [71 which presents gray-scale plots 
of the time-evolution for (a) the mass-loss rate and (b) 
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Fig. 5. — Gray-scale plots showing time evolution of mass-loss rate in units of the point-star CAK value. The corresponding rotation 
rate and final mass-loss rate are given above each plot. Time is in units of the flow time t = i?* /vx, . Note in each the eventual constancy 
of the mass-loss rate with both radius and time, indicating relaxation to a steady-state solution. Note that the ranges for both gray-scale 
and time differ for each panel. 
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TABLE 1 
Comparison of Kink 
Locations 





f*kink, hydro 

(R*) 


^ki nk. nozzle 

(«.) 


0.80 


4.18 


4.54 


0.82 


2.95 


3.03 


0.83 


2.52 


2.54 


0.84 


2.12 


2.16 


0.85 


1.77 


1.86 



Note. — n<mk,hydro and 
r kink nozzle are the kink loca- 
tions found respectively from nu- 
merical hydrodynamical models 
and from the semi-analytic noz- 
zle analysis. 

the scaled velocity gradient in this lo = 0.85 case. An 
overloading of the wind and eventual collapse of the time- 
dependent solution are shown in these figures. In the 
first 300 flow times, there is an initial kink formation, 
but at its lower radius of r^ink = 1.77ii* the kink outflow 
speed is quite low, i>kink ~ 370 km/s, well below the local 
escape speed t'esc^kink) ~ 640 km/s. With the reduced 
line-driving in the decelerating region, the wind outflow 
now stagnates at a finite radius, r stag ~ 10i?». There 
material accumulates until it is eventually pulled back 
by the stellar gravity into a reaccretion onto the star, 
effectively quenching both the kink and base outflow. 

For this particular case, the outflow never fully recov- 
ers from this quenching, but for only slightly more rapid 
rotation, lo > 0.86, any flow stagnation occurs relatively 
close to the star, with a correspondingly faster and less 
massive reaccretion, followed by a recovery to a slow ac- 
celeration solution, as indicated in figures [3] and UJ In the 
very rapid, near-critical rotation case to = 0.90, no kink 
or flow stagnation forms, and the solution relaxes more 
directly to the slow solution, as illustrated in the lower 
right panel of figure [5) 

Finally, table 1 compares the radius of the kinks found 
in our numerical simulations with the radius at which 
the nozzle analysis indicates a steep acceleration solution 
can no longer be maintained. For cases in which the 
mass flux set at the base rh — rho — n{x = 0) is above 
the CAK value, i.e. rh — n(x = 0) > 1 , this occurs 
at a radius Xk where a declining nozzle function falls 
back to n(xk) = rh. The corresponding radius rk; n k. nozzle 
agrees quite well with the kink locations found from the 
hydrodynamical simulations, r^ink, hydro- Note also that 
the requirement n(x = 0) — (4/9)/(l — lo 2 ) > 1 implies 
rotation rates of lo > V5/3 ~ 0.745, representing the 
onset for either possible kink solutions or a switch to a 
slow acceleration solution. 

5.3. Results for Slow Acceleration Initial Condition 

A central finding of the dynamical simulations is that 
moderately fast rotation models 0.75 < u> < 0.85 form 
fast-acceleration kink solutions instead of the slow accel- 
eration solutions predicted from steady-state analyses. 
But since the non-linear character of the flow equations 
allows more than one solution, this raises the question of 
whether slow acceleration solutions in this regime might 
also be stable attractors, perhaps for initial conditions 
that are closer to their slower outflow form than the fast, 



non-rotating model used for the initial conditions in the 
simulation models discussed above. For each of the spe- 
cific rotation cases in this transitional range, lo = 0.75, 
0.80, 0.82, and 0.84, we thus recompute simulations that 
instead use an initial condition set to the slow accelera- 
tion steady-state found for the faster, near-critical rota- 
tion case lo = 0.9. 

For the lo = 0.75 case, we find that the model again 
relaxes to a fast solution with an outer-wind kink at 
Hank ~ 10/2* . Steep acceleration solutions are also re- 
covered in all the slower rotation models as well. 

However, as shown in figure [H for the lo — 0.80, 0.82, 
and 0.84 cases the final states now do approach slow ac- 
celeration solutions, apart from a persistent peculiar up- 
ward kink near the outer boundary (i.e. for x > 0.9 in 
figure [5J . Such models also show persistent small-scale 
fluctuations with an amplitude of ca. 10% in the mass 
flux, apparently reflecting some difficulty for the numeri- 
cal solution to relax to a subcritical flow solution in these 
cases. The upward kink may stem from using a super- 
critical outflow boundary condition for this slow acceler- 
ation solution, which does not become supercritical until 
far from the star. 

Finally, even for this slow acceleration initial condition, 
the peculiar case w=0.85 still forms an overloaded con- 
dition with flow stagnation and reaccretion, after which 
it never fully recovers a steady outflow result. 

5.4. ID Results for Equatorial Density 

These results also allow us to identify the radial vari- 
ation of the relative density enhancement in the slow 
equatorial wind of a rotating star, compared to the non- 
rotating solution that applies to the polar wind. We are 
ignoring here gravity darkening and oblateness effects, as 
well as any "bi-stability" in the line-driving parameters 
between the polar and equatorial wind. From the analy- 
sis of §[31 the relative density enhancements are given by 
the ratios of the quantity rh/ ^Jw between the rotating 
and non-rotating models. 

The dashed curves in figure [5] show the spatial vari- 
ation of this density enhancement for rotating models 
with lo = 0.8, 0.86, 0.9, and 0.95. Note that the en- 
hancements are a few factors of ten, not insignificant, 
but not sufficient to reproduce the inferred densities 
of B[e] disks, which are factors of orde r 10 4 or more 
denser than a typical polar wind outfl ow ( Zickgra f et al.1 
119851 : iKraus fc Miroshnichenkol |2006[ ) . Inclusion of bi- 
stability effects could give about another factor of a 
few by increasing the equatorial mass loss around the 
cooler, gravity-darkened equ ator, thus yielding an over- 
all enhancement of or der 10 2 (|Lamers & Pauldrach 1991; 
iPelupessv et al.ll2000l ) . Significantly higher enhancement 
would require an unrealistically low a to increase further 
the equatorial mass flux, and/or assuming an equato- 
rial surface rotation within a sound speed of the critical 
(orbital) speed. In the latter case, minor disturbances 
(e.g. pulsations) in the stellar envelope or photosphere 
could instead eject material into an orbiting, Keplerian 
disk (jLee et al.1 119911 : lOwockil l2005h , obviating the need 
to invoke any central role for radiatively driven outflow 
solutions. 

6. SUMMARY AND CONCLUSIONS 
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Fig. 7. — (a) Gray-scale plot of the local mass-loss rate in units of the point-star CAK value for u) = 0.85, showing the collapse of the 
solution. Note the pile-up of material at a distance of r ~ 10i?* , and the eventual stagnation of the wind and reaccreation of material back 
onto the stellar surface, (b) Gray-scale plot of the scaled velocity gradient in cgs units for u> = 0.85. In the time interval between t = 150 
and 300 flow times, a kink in the velocity gradient is clearly visible, with the kink location at r^ink 1.77-R*. After t « 300, this kink 
solution becomes overloaded and unstable, eventually collapsing. 
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Fig. 8. — Final wind velocity law (flow speed over escape speed, 
v/Vesc, plotted vs. scaled inverse radius x) for lj = 0.8, 0.82, and 
0.84 simulations that used the u) = 0.9 final state (dashed line) 
as an initial condition. The solid lines now approximate the ex- 
pected slow acceleration solution, except for an unexplained uptick 
in velocity gradient near the right boundary, x > 0.9. 
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Fig. 9. — Density enhancement of slow wind solutions with 
rapid rotation rates, lo = 0.8, 0.86, 0.9, and 0.95, relative to a 
non-rotating wind with the same wind parameters. Dashed curves 
are for analytic results, while solid curves are those obtained from 
numerical simulations using VH-1. 



Using an analytic approach combined with numerical 
hydrodynamic simulations, we investigate the reasons for 
the switch from a steep to shallow acceleration in ID 
line-driven stellar wind models as stellar rotation rates 
are increased beyond a threshold value of co f« 0.745. The 
results indicate that the cause of this switch is the over- 
loading of the base mass-loss rate beyond the point-star 
CAK value. The latter represents the maximal allowed 
mass loss for which there can be a monotonically accel- 
erating flow speed throughout the entire wind. Further- 
more, the finite-disk correction-factor (fdcf) reduces the 
driving effectiveness near the stellar surface, and thus re- 
duces the maximal mass loss that can be initiated there. 
This reduction allows the outer wind to maintain a pos- 
itive acceleration even as other effects (e.g. centrifugal 
reduction in effective gravity near the surface) allow for 
an increase of the base mass loss from its fdcf value. 

This problem of wind overloading at la r ge ro tation 
rates was first noticed by iFriend fc Abbott) ([1986) , who 
found indeed that beyond some threshold rotation rate, 
the super-critical, steep acceleration solutions they were 
deriving could not be followed beyond some finite ra- 
dius. The contributions of M. Cure and collaborators 
have since shown that this termination can be avoided 
by switching to a shallow acceleration solution. As rota- 
tion increases beyond the threshold value of lo ss 0.745, 
the base mass-loss becomes greater than the point-star 
CAK value, and so the only globally accelerating solution 
possible is a shallow one with a subcritical outflow. 

However, an important lesson learned from our numeri- 
cal simulations is that the flow does not necessarily follow 
the solution with a globally monotonic acceleration. We 
see that the steady ID solutions for rotating winds actu- 
ally fall into four domains. First, for no or low rotation 
rates (lo < 0.745), the minimum of the nozzle function is 
less than unity and occurs at the stellar surface (x = 0). 
This allows the flow to transition to a super-critical out- 
flow directly from the static surface boundary and follow 
the steep acceleration solution with a mass-loss rate that 
is less than the point-star CAK value. 

Next, there exits a "gray zone" for rotation rates be- 
tween co = 0.745 and ~ 0.86 where two different solutions 
are possible. If we restrict ourselves to a strictly positive 
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acceleration, w' > 0, then as rotation increases beyond 
id ss 0.745 the base mass loss exceeds the point-star CAK 
value, and so the only globally accelerating solution pos- 
sible is one with a shallow acceleration and a mass-loss 
rate that saturates to the point-star CAK value. On the 
other hand, if we provide a backup scaling for negative 
accelerations and investigate cases where multiple local 
minima appear in the nozzle function, then overloaded 
situations where the square root term in equation (|19| 
becomes negative lead to a "kink" , and thus a negative 
acceleration or coasting solution. Moreover, for these 
new kink solutions, the multiple minima of the nozzle 
function now give rise to mass-loss rates that exceed the 
point-star CAK value. 

This "gray zone" of kink solutions can furthermore be 
divided into two domains. For rotation rates between 
uj w 0.745 and ~ 0.85, the overloading of the wind is 
not too severe. The coasting solution can thus still reach 
large distances with a finite speed, and so, these cases 
are able to approach a perfectly steady, time-independent 
state, despite the presence of the discontinuity in slope 
at the kink. On the other hand, for rotation rates of 
u> = 0.85, the wind becomes strongly overloaded and 
eventually stagnates at some finite radius. After enough 
material has piled up at this radius, it will fall back as 
a "reaccreation front" toward the star. If the bound- 
ary conditions are such that this material is allowed to 
fall through the lower boundary, then it effectively dis- 
appears from the model allowing the simulation to re- 
establish a wind outflow with a slow acceleration solution 
(as in the cu = 0.86 case). 

For rotation rates greater than u> « 0.86, the global 
minimum of the nozzle function is unity and occurs at 
large radii (x = 1). In order to satisfy the static surface 
boundary condition, the flow at all finite radii should be 
subcritical (but still supersonic) , and so the only solution 
possible is the shallow (— ) acceleration solution, with a 
mass-loss rate that has saturated to the point-star CAK 
value. 

It is also worth noting that since the non-linear na- 
ture of the flow equations allow more than one solution, 
slow acceleration solutions in the 0.75 < u> < 0.85 regime 
are also possible via a change in initial conditions. By 
initializing the wind using the slow acceleration steady- 
state found for the u> = 0.9 case, it is possible to achieve 
final states that do approach the expected slow acceler- 
ation solution (apart from a peculiar upward kink near 
the outer boundary, see fig. UJ). 

One should also add here a few words regarding the 
validity of the Sobolev approximation for the slow accel- 
eration solutions. Recall that for the usual case where 
the broadening of the line is set by the local ion ther- 
mal speed Vth, the geometric width of the local res- 
onance with the r adiation is ab out a Sobolev length, 
Isob = v t h/{dv/dr) (|Sobolevill960l ). Thus, in a supersonic 
flow, this Sobolev length is of order v t h/v <C 1 smaller 
than a typical flow variation scale, such as the den- 
sity/velocity scale height H = \p/(dp/dr)\ w v/{dv/dr). 
Therefore, while the flow speeds for the slow acceleration 
solutions are lower than those of the steep acceleration 
solutions, they are still supersonic throughout nearly the 
entire wind, and so the Sobolev approximation remains 
quite appropriate. 



As noted in the introduction, it is important to keep 
in mind the limited physical relev ance of such ID mod - 
els. The 2D WCD si mulat ions of 
iCranmer fe Owockl ([1995), and 
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serve as an example, showing that for moderately rapid 
rotation, there can be a 2D flow pattern by which mate- 
rial from higher latitudes is focussed toward the equa- 
tor through the WCD effect. Depending on whether 
the material reaches the equator above or below some 
"stagnation point" , it either drifts outward or falls back 
toward the star. This simultaneous infall plus outflow 
behavior is not possible in a steady ID model, but 
is a perfectly natural occurrence in a 2D simulation. 
There are also the issues of oblateness, limb darkening, 
and the non-radial line force, and how these affect the 
latitudinal motion of the flow leading to inhibition of 
the WCD (ICranmer fc OwockiH995t lOwocki et alJll996t 
iPetrenz fc Pulsll2000f) . Finally, equatorial gravity dark- 
ening can reduce the wind mass flux from the equator 
and lead to an e quatorial wind densi ty that is lower than 
near the poles (|Owocki et al.l [19961) . A key point here 
is that ID simulations represent a sort of "best-case" 
scenario for the formation of a disk. A move to 2D sim- 
ulations shows that when the relevant physics (such as 
nonzero, nonradial line forces and gravity darkening) are 
included, material will tend to be channeled away from 
the equator and inhibit disk formation. Thus, if ID sim- 
ulations are incapable of producing equatorial densities 
capable of explaining those inferred in B[e] supergiants, 
then 2D simulations will most likely not change this. 

Beyond even developing a 2D CAK-type model, there 
remain key physical limitations not accounted for in this 
CAK formalism. Two examples are the intrinsic, small- 
scale instability of line-driving, and multiple scattering 
effects. The former might well disrupt a slow-acceleration 
solution, even though the CAK form of such solutions 
seem from the present simulations to be a stable attrac- 
tor. The latter places strict upper limits on the mass flux 
that can be driven within a geometrically thin disk, since 
even with ideally tuned line-driving parameters within a 
CAK model (e.g. choosing an anomalously small CAK 
exponent a to enhance the expected CAK-type mass loss 
rate), radiative driving in the thin disk would generally 
be limited to one or two scat terings before the photon es- 
capes out of the disk plane (|Owockill2006h . Thus, quite 
generally, radiative driving could not produce a disk out- 
flow that exceeds the single scattering limit. So, to the 
extent that observational inferences of supergiant B[e] 
disks imply a very dense medium (dense enough to form 
dust), it seems more likely that these represent orbiting 
Keplerian disks, with perhaps ablation flows off the disk 
surface producing the equatorial outflow inferred by ob- 
servations of Doppler-shifted absorption troughs in UV 
resonance lines. 
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APPENDIX 

OUTER CRITICAL POINT LOCATION FOR SHALLOW SOLUTIONS 
WITH FINITE SOUND SPEED 

Let us examine here a key issue regarding the location of the outer critical point for the shallow acceleration solutions 
that appear for rotation rates lu > 0.75. For the zero sound-speed limit (as analyzed in §§2-3), the radius of this critical 
point formally approaches infinity, r c — > oo, or x c — > 1 in the scaled inverse radius coordinate x = 1 — i?»/r. This 
raises the issue of whether the numerical hydrodynamic simulations (§§4-5) done with a finite outer boundary radius 
must use some sort of subcritical right boundary condition (with one inward and one outward pointing characteristic), 
instead of the supercritical outflow boundary conditions (with two outward characteristics) that are usually used in 
stellar wind simulations. 

Note however that the hydrodynamic simulations presented above include a small but non-zero sound speed, char- 
acterized by the dimcnsionless parameter w s = a 2 /v^ sc i=s 4 x 10~ 4 . When this is included in the equation of motion 
(as in eq. [5]), the additional term on the right side has the scaled form 4w s /(l — x) (we are neglecting the additional 
w s /w term on the left hand side of eq. [6]). The total effective gravity used in the denominator of the nozzle function 
(see eq. [20]) then takes the form 

g (x) = l-u\l-x)-^-. (Al) 
1 — x 

Instead of declining monotonically as x increases, this now has a maximum at some finite radius, and that in turn 
allows the nozzle function to reach a minimum at some finite radius. In a model with a non-zero sound speed, this 
thus becomes the effective critical point of the slow acceleration solutions. Since w s -C 1, this minimum still occurs 
near x < 1, and so to evaluate the location, it is convenient to expand the finite disk correction factor about x = 1. 
For a j3 = 1 velocity law and a = 1/2, this gives 

; s-l, (A2) 

which makes the nozzle function 

nix) = Z!M ~ l + (l-s)/2 

1 ' ~ g(x) 1 - -x)- 4 Ws /(l - x) ■ [ 6) 

We can easily solve for the critical location x c through 

n'(x c ) = Q. (A4) 

To lowest order in the small parameter y/w^ ~ 0.02, we find 

So «l- 2 y=^:«l--^ 7f «l-0.04VS(l- W a ), (A5) 
which implies a critical radius of 



R* 



17.67 (1 + w 2 ). (A6) 



The upshot of this is thus that the outer boundary radius used in our numerical models (i? max = 100i?„) should 
be well above the critical point (r c < 36i?*) for even the shallow slope solutions at high rotation rate. Thus a pure 
super-critical outflow boundary condition should indeed be appropriate for all the models computed here. 
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